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Abstract 



["t"] ■ This paper develops a bias correction scheme for a multivariate heteroske- 

dastic errors-in-variables model. The applicability of this model is justified in 
areas such as astrophysics, epidemiology and analytical chemistry, where the 
variables are subject to measurement errors and the variances vary with the 
observations. We conduct Monte Carlo simulations to investigate the perfor- 
mance of the corrected estimators. The numerical results show that the bias 
correction scheme yields nearly unbiased estimates. We also give an applica- 
tion to a real data set. 
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1 Introduction 



Heteroskedastic errors-in-variables (or measurement error) models have been exten- 
sively studied in the statistical literature and widely applied in astrophysics (to 
explain relationships between black hole masses and some variates of luminosities), 
epidemiology (to model the cardiovascular event with its risk factors), analytical 
chemistry (to compare different types of measurement instruments). The applica- 
bility of this model abound mainly in the astronomy literature where all quantities 
are subject to measurement errors (Akritas and Bershady, 1996). 

It is well-known that, when the measurement errors are ignored in the estimation 
process, the maximum-likelihood estimators (MLEs) become inconsistent. More 
precisely, the estimation of the slope parameters is attenuated (Fuller, 1987). When 
variables are subject to measurement errors, a special inference treatment must 
be carried out for the model parameters in order to avoid inconsistent estimators. 
Usually, a measurement equation is added to the model to capture the measurement 
error effect and then the MLEs from this approach are consistent, efficient and 
asymptotically normally distributed. A careful and deep exposition on the inferential 
process in errors-in-variables models can be seen in Fuller (1987) and the references 
therein. 
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Although consistent, asymptotically efficient and asymptotically normally dis- 
tributed, the MLEs are oftentimes biased and point inference can be damaged. This 
is not a serious problem for relatively large sample sizes, since bias is typically of 
order (9(n _1 ), while the asymptotic standard errors are of order (9(n -1 / 2 ). How- 
ever, for small or even moderate values of the sample size n, bias can constitute a 
problem. Bias adjustment has been extensively studied in the statistical literature. 
For example, Cook et al. (1986), Cordeiro (1993), Cordeiro and Vasconcellos (1997), 
Vasconcellos and Cordeiro (1997) and, more recently, Cordeiro (2008). Additionally, 
Patriota and Lemonte (2009) obtained general matrix formulae for the second-order 
biases of the maximum-likelihood estimators in a very general model which includes 
all previous works aforementioned. The model presented by the authors considers 
that the mean vector and the variance-covariance matrix of the observed variable 
have parameters in common. This approach includes the heteroskedastic measure- 
ment error model that we are going to study in this paper. 

The main goal of this article is define bias-corrected estimators using the general 
second-order bias expression derived in Patriota and Lemonte (2009) assuming that 
model defined by (1) and (2) holds. Additionally, we compare the performance of 
bias-corrected estimators with the MLEs in small samples via Monte Carlo simu- 
lations. The numerical results show that the bias correction is effective in small 
samples and leads to estimates that are nearly unbiased and display superior finite- 
sample behavior. 

The rest of the paper is as follows. Section 2 presents the multivariate het- 
eroskedastic errors-in-variables model. Using general results from Patriota and 
Lemonte (2009), we derive in Section 3 the second-order biases of the MLEs of 
the parameters. The result is used to define bias-corrected estimates. In Section 4 
the 0(n- v ) bi ases of the estimates fa and S« are given. Monte Carlo simulation 
results are presented and discussed in Section 5. Section 6 gives an application. 
Finally, concluding remarks are offered in Section 7. 



2 The model 

The multivariate model assumed throughout this paper is 

Vi = fa + fcXi + q u i = l,...,n, (1) 

where t/i is a (v x 1) latent response vector, Xi is a (m x 1) latent vector of covariates, 
(3 is a (v x 1) vector of intercepts, fli is a (v x m) matrix which elements are 
inclinations and is the equation error having multivariate normal distribution 
with mean zero and covariance- variance matrix S g . The variables yi and Xi are 
not directly observed, instead surrogate variables Yi and Xi are measured with the 
following additive structure: 

Y i = y i + riy i and X i = x i + ri a> .. (2) 

The errors r\ yi and t] Xi are assumed to follow a normal distribution given by 

te) ~M(o) ■ ft :)]■ 
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where "~" means "independently distributed as" and the covariance- variance ma- 
trices T yi and T Xi are assumed to be known for alH = 1, . . . , n. These matrices may 
be attained, for example, through an analytical treatment of the data collection 
mechanism, replications, machine precision, etc. 

Model (2) has equation errors for all lines, i.e., %/i and Xi are not perfectly related. 
These equation errors are justified by the influence of other factors than Xi in the 
variation of y L . It is very reasonable to consider equation errors in (1) to capture 
extra variability, since the variances r Vi are fixed and whether some other factor 
affects the variation of the estimation of the line parameters will be clearly 
affected. Supposing that x^ ~ Af m (n x , SJ, where "~" means "independent and 
identically distributed as", and considering that the model errors (qi, r\ yi and rj Xi ) 
and Xi are independent, we have that the joint distribution of the observed variables 
can be expressed as 

(I) 

Note that in (3), the mean vector and the covariance- variance matrix of observed 
variables have the matrix (3\ in common, i.e., they share mv parameters. Kulathinal 
et al. (2002) study the univariate case (when v — 1 and m — 1) and propose an EM 
(Expectation and Maximization) algorithm to obtain MLEs for model parameters. 
In addition, they derived the asymptotic variance of the MLE of the inclination 
parameter making it possible to build hypotheses testing of it. Also, de Castro et 
al. (2008) derive the observed and expected Fisher information and conduct some 
simulation studies to investigate the behavior of the likelihood ratio, score, Wald and 
C(pt) statistics for testing hypothesis of the parameters and Patriota et al. (2009) 
study the asymptotic properties of method-of-moments estimators in the univariate 
model proposed by Kulathinal et al. (2002). Model (2) is a multivariate version of 
the model proposed by Kulathinal et al. (2002). 

3 Second-order bias of 6 

In order to follow the same scheme adopted by Patriota and Lemonte (2009), define 
the vector of parameters = (0J , vec(/3i) T , /x^, vech(£,j,) T , vech(S 9 ) T ) T , where 
vec(-) is the vec operator, which transforms a matrix into a vector by stacking 
the columns of the matrix and vech(-) is the vech operator, which transforms a 
symmetric matrix into a vector by stacking the on or above diagonal elements. 
Also, consider Zj = (Yj , Xj) T and the mean and covariance-variance function as 

* W =(A + AM and ^e)=(^i + J + T " *M, 

respectively. 

Moreover, to simplify notation, define the quantities Z = vec(Zi, . . . , Z n ), fi = 
vec(Mi(0),...,Mn(0)), £ = diag{S 1 (6>),...,S n (6>)} and u = Z - fi. The log- 
likelihood function for the vector parameter 6 from a random sample, except for 



A) + 01 Ma 
^x 



P^xPj + S g + r m (3{E a 



Sec + T x 



(3) 
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constants, can be expressed as 

= — log |E|-itr{E-W}. (4) 

Additionally, for the purpose of computing the score function, the Fisher information 
and the second-order biases, also define 

d» gV C - &Cr A - -Y,- l C E- 1 

r ~ dd r : sr ~dd s ddr- de r ' sr ~de s ' r ~ r 



and 



df3 , a) dp! _ (a) duo 



(r) = U^JO (s) = UJJi (s) _ _ 

„( s ) dT, x ( s ) <9E 9 

F *- = «r and < 

with r, s = l,2,...,p, where p is the dimension of 0. The quantities (5) are 
vectors or matrices of zeros with a unit in the position referring to the s th ele- 
ment of 0. Let D = (0^,0^,0^, 0,0) and V = (0, C^, 0, C Sa! , C S J, with 

a A) = ( a l) • • • ) a v)i = ( a D+l) • • • ) a v(m+l)), a ^ — ( a v(rn+l)+l, ■ ■ ■ , a v(m+l)+m), 

C/3i = fvec(C w+ i),...,vec(C w ( m+ i))), C Sa; = (vec(C ( „ + i)( m+ i)), . . . , vec(Cy)) and 
Cs q = (vec(Cy+i), . . . , vec(Cp)), where p' = v(m + 1) + m + m(m + l)/2. 
The first derivative of (4) with respect to the r th element of 6 is 

U r = itr{A r (E - uu T )} + tr{E - Vw T }; (6) 

the expectation of the derivative of (6) with respect to the s th element of 6 is given 
by 

n sr = -tr{A r C s } — ajE -1 a r . 
2 

Note that, under general regularity conditions (Cox and Hinkley, 1974, Ch. 9), 
— K sr is the (s,r) th element of the expected Fisher information. The score function 
and the expected Fisher information are given, respectively, by Ug = D T Y,~ 1 u — 
|V T E- 1 vec(E - uu T ) and K e = D T T,~ 1 D + |y T E" 1 V, with E = E <g> E and <g> 
is the Kronecker product. Defining 

*=( ,v" A and H=(l X" 

V-vec(E - uu T )J ' Ivy V° 2E 

we can write the score function and the Fisher information in a short form as 

U = FHu and K e = FHF. 

The Fisher scoring method can be used to estimate iteratively solving the 
equation 

0(m+l) = ^F(rn)Tjj(m)p(m)ylp(m)Tjj(m)~*(m)^ m = 0,1,2, ... , (7) 
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where u*( m ^ = F( m ) Q( m ) -\-u( m \ Each loop, through the iterative scheme (7), consists 
of an iterative re- weighted least squares algorithm to optimize the log-likelihood (4). 
Using equation (7) and any software (MAPLE, MATLAB, Ox, R, SAS) with a weighted 
linear regression routine one can compute the MLE, 0, iteratively. Initial approx- 
imation for the iterative algorithm is used to evaluate F^°\ and u*^ 
from which these equations can be used to obtain the next estimate 0^ l \ This new 
value can update F, H and u* and so the iterations continue until convergence is 
achieved. 

The general matrix formulae derived by Patriota and Lemonte (2009) for n~ l 
bias vector B(0) of is given by 



B{9) = (F T HF)- 1 F T H^, 



(8) 



where £ = . . . , $ p )vec{(F T HF)' 1 } and $ r = -\{G r + J r ), r = 1, 2, 



with 



G r 



vec(Ci r ) 



a 



pr 



vec(C, 



pr j 



and J r 







2(1 



nq 



a r )D 



where Ik denotes an k x k identity matrix. The bias vector B(0) is simply the 
set coefficients from the ordinaryjweighted lest-squares regression of the £ on the 
columns of F, using weights in H . The bias vector B(0) will be small when £ 
is orthogonal to the columns of H and it can be large when n is small. Note that 
equation (8) involves simple operations on matrices and vectors and we can calculate 
the bias B(0) numerically via software with numerical linear algebra facilities such 
as Ox (Doornik, 2006) and R (R Development Core Team, 2006) with minimal effort. 
After a somewhat algebra, we have 



A) 









Of 



? (t) 



and a, 



0, 



for r = 1, . . . , v; s — v + 1, . . . , v(m + 1); t — v(m + 1) + 1, . . . , v(m + 1) + m; and 
u — (v + l)(m + 1), . . . ,p; where p = v(m + 1) + m + m(m + l)/2 + v(y + l)/2. 
(Here, l n denotes an n x 1 vector of ones.) Moreover, 



a r 







for all r and s, 



c s -i n ®^ F ^ q j,C t -I n 



and 



C — T 



F< u) 



0/ ' 



F%iPl 
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for s — v + 1, . . . , v(m + 1); t — v(m + 1) + 1, . . . , v(m + 1) +m; and u — (v + l)(m + 
l),...,p. Additionally, 

C rs = I n <g> [ * A x pi + Pi x Pi 







and 



c -J J^i+^K" 



S, 1 /3i 

for r, s, -u = v + 1, . . . , i>(m + 1); £ = v(m + 1) + 1, . . . , v(m + 1) + to; and C rs = 
otherwise. 

Therefore, in the measurement error model defined by the equations (1) and 
(2), all quantities necessary to compute the 0(n _1 ) bias of 6 using expression (8) 
are given. On the right-hand side of expression (8), consistent estimates of the 
parameter 6 can be inserted to define the corrected MLE 6 = 6 — B(0), where B(-) 
denotes the MLE of B(-), that is, the unknown parameters are replaced by their 
MLEs. The bias-corrected estimate 6 is expected to have better sampling properties 
than the uncorrected estimator, 6. In fact, we present some simulations in Section 5 
to show that 6 has smaller bias than its corresponding MLE, thus suggesting that 
the bias corrections have the effect of shrinking the modified estimates toward to 
the true parameter values. 



4 Biases of the MLEs jli and Ej 

In this section, we give matrix formulae for the 0(n _1 ) biases of the MLEs of the 
ith mean \Xi = fXi(6) and ith variance-covariance vector X* = vech(Ej(0)). Let 
q± — v + m and qi = q\{q\ + l)/2. Additionally, let A = [A±, . . . , A n ] T be a np x p 
matrix, where A; is a p x p matrix, then we define tr*(A) = [tr(Ai), . . . , tr(A n )] T . 

From a Taylor series expansion of fii = /J>i(6), we obtain up to an error of order 
0(n- 2 ): 

B(fii) = LiBfi) + hr*[M t Cov(6)], 

where Li is a q\ x p matrix of first partial derivatives d\XijdQ r (for r — 1, 2, . . . 
Mi = [Mn, . . . , M iqi ] T is a q±p x p matrix of second partial derivatives, where Mu 
is a p x p matrix with elements d 2 fiu/d9 r d9 s (for r, s — 1, . . . ,p and I = 1, 2, . . . , qi), 
Cov(0) = Kq 1 is the asymptotic covariance matrix of 6 and the vector B{6) was 
defined before. All quantities in the above equation should be evaluated at 6. The 
asymptotic variance of /Xj can also be expressed explicitly in terms of the covariance 
of 6 by 

Var(^) = LiCov(6)Lj. 

The second-order bias of S* is obtained by expanding S* = £*(#) in Taylor 
series. Then, the 0(n _1 ) bias of S* is written as: 

B(S*) = L*B(6) + itr*[M*Cov(0)], 
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where L* is a q 2 x p matrix of first partial derivatives d"S*/d6 r (for r = 1, 2, . . . ,p), 
M* = [M*!, . . . , M* ] T is a g 2 £> x p matrix of second partial derivatives, where M.* d 
is a p x p matrix with elements 5 2 E* Z / 'd9 r d9 s (for r, s = 1, . . . ,p and / = 1, 2, . . . , q 2 ). 

Therefore, we are now able to define the following second-order bias-corrected 
estimators for ft and £*: 

ft = ft - B(ft) and E* = E*-B(£*). 

It is clear that the 0(n~ l ) bias of any other function of 6, say ^(6) (h x 1), can 
be obtained easily by Taylor series expansion: 

= vg>B(0) + VfV^Cov^)], 

where is a h x p matrix of first partial derivatives d^/d6 r (for r = 1, 2, . . . ,p) 
and = [V^, . . . , V^)j T is a hp x p matrix of second partial derivatives, where 
V^) is a p x p matrix with elements d 2 ^i/d9 r d9 s (for r, s = l,...,p and / = 
1,2,..., h). 



5 Numerical results 

We shall use Monte Carlo simulation to evaluate the finite sample performance of 
the MLEs attained using the iterative formula (7) and of their correspondent bias- 
corrected versions for a heteroskedastic errors-in-variables model presented in (2) 
with m — v — 1. The sample sizes considered were n = 40,60,100 and 200, the 
number of Monte Carlo replications was 10,000. All simulations were performed 
using the R programming language (R Development Core Team, 2006). 
We consider the simple errors-in-variables model 

Yi = Hi + Vy t and Xi = Xi + r) Xi , 

with yi\xi ~ d A/"(A) + Pi^i, <y 2 ). This model was studied by Kulathinal et al. (2002). 
The errors r] y . and r\ x . are independent of the unobservable covariate Xi and are 
distributed as 

W N 





ind 


w 




[ r T . 


are 



M 2 



that Xi ~ J\f(n x ,al), we have that the joint distribution of the observed variables 
can be expressed as 











x% 


1 % 


= 1, 



i \ ind .r 



+ r yi + a 2 Pi<J x 



< + Tx 



Define 6 = (#,, /3 U fi x , a 2 x , a 2 ) J , 

» l (0)=( (30 + (3lfix ) and £,(0) 
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From the previous expressions, we have immediately that 



ai = l n (8) I 1 J , a 2 = l n ® J , a 3 = l n <g>(f), a 4 = a 5 = 



v oy ' z n ~ v o j ' a n ~ v i 

and a rs = for all r, s except for 



a 23 — a 32 — In 



Also, d = C 3 = and 



C* 2 = /„ ( 2P f» i),C 4 = I n ®( P * f| and /„ ' ' ° 



Additionally, C rs = for all r, s except for 

C 22 = I n ® ^ ^ and C 24 = C 42 = I„ ® ^ J 

Thus, I? = (a 1; a 2 , a 3 , 0, 0) and V = (0, vec(C 2 ), 0, vec(C 4 ), vec(C 5 )). Therefore, 
all the quantities necessary to calculate B(0) using expression (8) are given. 

In order to analyze the point estimation results, we computed, for each sample 
size and for each estimator: relative bias (the relative bias of an estimator 9 is defined 
as [E(9) — 9}/9, its estimate being obtained by estimating E(0) by Monte Carlo) and 
root mean square error, i.e., VMSE, where MSE is the mean squared error estimated 
from the 10,000 Monte Carlo replications. For practical reasons and without loss 
of generality, we adopt the same setting of parameters chosen by de Castro et al. 
(2008). (The parameters are the MLEs for the model parameters using a real data 
set presented in the next section.) We take (3 = —2, j3i = 0.5, /i x = —2, a 2 £ = 4 and 
a 2 = 10. We also consider two types of heteroskedasticity as studied by Patriota 
et al. (2009), namely: (a) v /t^ ~ [7(0.5,1.5) and jr^ i ~ [7(0.5, 4), where U(a,b) 
means uniform distribution on [a,b]; (b) ^J~f^ i = 0.1\xi\ and y/Ty~ = 0.1| — 2 + 0.51Xj|, 
i.e., the variances depend on the unknown covariate. We remark that the variances 
are considered to be known and keep fixed in all Monte Carlo simulations. 

Table 1 shows simulation results for an errors-in-variables model with a uniform 
heteroskedasticity. The figures in this table reveal that the maximum-likelihood 
estimators of the parameters can be substantially biased when the sample size is 
small, and that the bias correction we derived in the previous section is very effective. 
For instance, when n = 40 the biases of the estimators of /3q, /3\, fi x , cr 2 and a 2 
average —0.02244 whereas the biases of the corresponding bias-adjusted estimators 
average —0.00276; that is, the average bias (in value absolute) of the MLEs is almost 
ten times greater than that of the corrected estimators. In particular, the maximum- 
likelihood estimators of a 2 and a 2 display substantial bias, and the bias correction 
proves to be quite effective when applied to these estimators. 

Table 2 displays simulation results for an errors-in-variables model with a nonuni- 
form heteroskedasticity. We note that the bias-adjusted estimator again displays 
smaller bias than the standard maximum-likelihood estimator. This suggests that 
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the second-order bias of MLEs should not be ignored in samples of small to mod- 
erate sizes since they can be nonnegligible. Note also that root mean square error 
decrease with n, as expected. Additionally, we note that all estimators have similar 
root mean squared errors. 

It is interesting to note that the finite-sample performance of the estimator of 
a 2 deteriorate when we pass of the model with a uniform heteroskedasticity for the 
model with a nonuniform heteroskedasticity (see Tables 1 and 2). For instance, 
when n = 100, the relative biases of a 2 (MLE) were —0.0135 (uniform heteroskedas- 
ticity) and —0.0484 (nonuniform heteroskedasticity), which amounts to an increase 
in relative biases of nearly 3.5 times. 

Table 1: Relative bias and VMSE of uncorrected and corrected estimates with a 
uniform heteroskedasticity: ^/f^ ~ [7(0.5, 1.5) and ^fr^ ~ £7(0.5,4). 

MLE BCE 



n 





Rel. bias 


VMSE 


Rel. bias 


VMSE 


40 


Po 


-0.0173 


0.99 


-0.0043 


0.97 




Pi 


0.0315 


0.38 


0.0054 


0.37 






-0.0018 


0.35 


-0.0018 


0.35 




°l 


-0.0351 


1.11 


-0.0045 


1.13 




a 2 


-0.0895 


3.31 


-0.0086 


3.38 


60 


Po 


-0.0139 


0.77 


-0.0061 


0.76 




Pi 


0.0213 


0.29 


0.0058 


0.29 




Hx 


0.0009 


0.28 


0.0009 


0.28 




-I 


-0.0239 


0.89 


-0.0036 


0.90 




a 2 


-0.0548 


2.60 


-0.0018 


2.64 


100 


Po 


-0.0100 


0.68 


-0.0037 


0.67 




Pi 


0.0168 


0.26 


0.0042 


0.25 






0.0001 


0.25 


0.0001 


0.25 




°l 


-0.0135 


0.80 


0.0022 


0.81 




a 2 


-0.0424 


2.40 


0.0003 


2.43 


200 


Po 


-0.0049 


0.59 


-0.0006 


0.59 




Pi 


0.0127 


0.22 


0.0041 


0.22 




Hx 


0.0013 


0.23 


0.0013 


0.23 




°l 


-0.0116 


0.70 


0.0008 


0.70 




a 2 


-0.0350 


2.09 


-0.0014 


2.11 



BCE: bias-corrected estimator. 
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Table 2: Relative bias and vMSE of uncorrected and corrected estimates with a 
nonuniform heteroskedasticity: = 0.1|xj| and ^Jr^ = 0.1\(3 + P\Xi\. 



MLE BCE 



Ti 


u 


rvei. uias 


V iv±oi-j 


rvei. uias 


V iv±oj_j 


a n 


a 

Po 


— U.UUzo 


U. ( O 


— U.UUlo 


n 7Q 
v. 1 6 




Pi 


n nooo 
u.uzyz 


U.Z I 


U.Uz I 


n 07 

U.Z ( 






n noos 

— U.UZZo 


n qo 


n noos 
— U.Uzzo 


n qo 




_2 
°x 


— U.Uoy4 


n oi 


— U.UO04 


n oo 




a 2 


-0.0540 


2.26 


-0.0056 


2.30 


60 


Po 


0.0008 


0.59 


0.0013 


0.59 




Pi 


0.0203 


0.22 


0.0192 


0.22 






-0.0208 


0.26 


-0.0208 


0.26 






-0.0502 


0.76 


-0.0340 


0.75 




a 2 


-0.0332 


1.85 


-0.0002 


1.88 


100 


Po 


0.0013 


0.51 


0.0016 


0.51 




Pi 


0.0184 


0.19 


0.0176 


0.19 






-0.0198 


0.23 


-0.0198 


0.23 




°l 


-0.0484 


0.65 


-0.0363 


0.65 




2 

a 


-0.0223 


1.61 


0.0027 


1.64 


200 


Po 


0.0036 


0.45 


0.0039 


0.45 




Pi 


0.0165 


0.17 


0.0159 


0.17 




Hx 


-0.0186 


0.20 


-0.0186 


0.20 






-0.0474 


0.59 


-0.0377 


0.58 




a 2 


-0.0204 


1.41 


-0.0004 


1.43 



BCE: bias-corrected estimator. 
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6 Application 



We shall now present an application of the model described in Section 2 where 
v — m — 1. We analyze a epidemiological data set from the WHO MONICA 
(World Health Organization Multinational MONitoring of trends and determinants 
in CArdiovascular disease) Project. This data set was previously studied by Kulathi- 
nal et al. (2002) and de Castro et al. (2008) where the ML approach was adopted 
to estimate the model parameters. 

The main goal of this project is to monitor trends in cardiovascular diseases (y) 
and relate it with known risk factors (x). The latent variables y and x are linearly 
related as 

Vi = A + PiXi + ft, i = 1, ■ ■ ■ , n. 

As the variables and Xi are not directly observable, surrogate variables and Xj 
are observed in their place, respectively. Such surrogate variables are attained from 
an analytical treatment of the data collection process. The data set are divided into 
two groups, namely: men (n = 38) and women (n = 36). 

In what follows, we compare the MLEs with the bias-corrected estimators. Ta- 
ble 3 presents the MLEs, its standard deviation, its second-order biases and the 
corrected estimates. It can be seen that, the greater is the standard deviation of 
the MLE, the more distant from zero is its respectively second-order bias. As con- 
cluded in the simulation studies, the biases of the variances estimates are larger 
than of those produced by the line estimators. The second-order biases of the MLEs 
can be expressed as a percentage of the MLEs. That is, for the men data set, the 
second-order biases are -0.21%, 0.85%, 0.00%, -2.92% and -9.21% of the total 
amount of the MLEs of f3 , A> fi x , a 2 and a 2 , respectively. For the women data 
set, the second-order biases are 52.96%, 1.21%, 0.00%, -3.16% and -10.19% of the 
MLEs of /3i, fi x , a 2 and a 2 , respectively. It shows that the second-order biases 
of the MLEs are more pronounced in the women data set, mainly for the intercept 
estimator. 



Table 3: MLEs and bias-corrected estimates. 





Parameter 


MLEs 


S.E. 


Bias 


BCEs 




A 


-2.0799 


0.5285 


0.0044 


-2.0843 




A 


0.4690 


0.2339 


0.0040 


0.4650 


Men 




-1.0924 


0.3550 


0.0000 


-1.0924 






4.3163 


1.0969 


-0.1261 


4.4423 




a 2 


4.8883 


1.7790 


-0.4501 


5.3384 




Parameter 


MLEs 


S.E. 


Bias 


BCEs 




A 


0.0321 


1.1121 


0.0170 


0.0151 




A 


0.6790 


0.4072 


0.0082 


0.6708 


Women 




-2.0677 


0.3386 


0.0000 


-2.0677 






3.6243 


0.9695 


-0.1146 


3.7389 




a 2 


11.0809 


4.2425 


-1.1289 


12.2098 



BCE: bias-corrected estimates. 
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7 Conclusions 



We derive a bias-adjustment scheme to eliminate the second-order biases of the 
maximum-likelihood estimates in a heteroskedastic multivariate errors-in-variables 
regression model using the general matrix formulae for the second-order bias derived 
by Patriota and Lemonte (2009). The simulation results presented show that the 
MLEs can be considerably biased. The bias correction derived in this paper is very 
effective, even when the sample size is large. Indeed, the bias correction mechanism 
adopted yields modified maximum-likelihood estimates which are nearly unbiased. 
Additionally, many errors-in-variables models are special cases of the proposed model 
and the results obtained here can be easily particularized to these submodels. We 
also present an application to a real data set. 
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